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Abstract. We consider the class of nonlinear optimal control problems (OCP) 
with polynomial data, i.e., the differential equation, state and control con- 
straints and cost are all described by polynomials, and more generally for 
OCPs with smooth data. In addition, state constraints as well as state and/or 
action constraints are allowed. We provide a simple hierarchy of LMI (lin- 
ear matrix inequality)-relaxations whose optimal values form a nondecreasing 
sequence of lower bounds on the optimal value. Under some convexity assump- 
tions, the sequence converges to the optimal value of the OCP. Preliminary 
results show that good approximations are obtained with few moments. 



1. INTRODUCTION 

Solving a general nonlinear optimal control problem (OCP) is a difficult chal- 
lenge, despite powerful theoretical tools are available, e.g. the maximum principle 
and Hamilton- Jacobi-Bellman (HJB) optimality equation. The problem is even 
more difficult in the presence of state and/or control constraints. State constraints 
are particularly difficult to handle, and the interested reader is referred to Capuzzo- 
Dolcetta and Lions [7] and Soner [33] for a detailed account of HJB theory in the 
case of state constraints. There exist many numerical methods to compute the so- 
lution of a given optimal control problem; for instance, multiple shooting techniques 
which solve two-point boundary value problems as described, e.g., in |40[ I34j . or 
direct methods, as, e.g., in [41 [ I12 | H5]. which use, among others, descent or gradient- 
like algorithms. To deal with optimal control problems with state constraints, some 
adapted versions of the maximum principle have been developed (see [251 [32] , and 
see [14] for a survey of this theory), but happen to be very hard to implement in 
general. 

On the other hand, the OCP can be written as an infinite-dimensional linear 
program (LP) over two spaces of measures. This is called the weak formulation of 
the OCP in Vinter [44] (stated in the more general context of differential inclusions). 
The two unknown measures are the state-action occupation measure (o.m.) up to 
the final time T, and the state o.m. at time T. The optimal value of the resulting 
LP always provides a lower bound on the optimal value of the OCP, and under 
some convexity assumptions, both values coincide; see Vinter |44j and Hernandez- 
Hernandez et al. 231 as well. 
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The dual of the original infinite dimensional LP has an interpretation in terms of 
" subsolutions" of related HJB-like optimality conditions, as for the unconstrained 
case. The only difference with the unconstrained case is the underlying function 
space involved, which directly incorporate the state constraints. Namely, the func- 
tions are only defined on the state constraint set . 

An interesting feature of this LP approach with o.m.'s is that state constraints, 
as well as state and/or action constraints are all easy to handle; indeed they simply 
translate into constraints on the supports of the unknown o.m.'s. It thus provides 
an alternative to the use of maximum principles with state constraints. However, 
although this LP approach is valid for any OCP, solving the corresponding (infinite- 
dimensional) LP is difficult in general; however, general LP approximation schemes 
based on grids have been proposed in e.g. Hernandez and Lasserre |21j . 

This LP approach has also been used in the context of discrete-time Markov 
control processes, and is dual to Bellman's optimality principle. For more details the 
interested reader is referred to Borkar [4 , Hernandez-Lerma and Lasserre [T8lH~9l[22] 
and many references therein. For some continuous-time stochastic control problems 
(e.g., modeled by diffusions) and optimal stopping problems, the LP approach has 
also been used with success to prove existence of stationary optimal policies; see for 
instance Cho and Stockbridge [5], Helmes and Stockbridge [15] , Helmes et al. |16j . 
Kurtz and Stockbridge [17] > an d also Fleming and Vermes [TT]. In some of these 
works, the moment approach is also used to approximate the resulting infinite- 
dimensional LP. 

Contribution. In this paper, we consider the particular class of nonlinear 
OCP's with state and/or control constraints, for which all data describing the prob- 
lem (dynamics, state and control constraints) are polynomials. The approach also 
extends to the case of problems with smooth data and compact sets, because poly- 
nomials are dense in the space of functions considered; this point of view is detailed 
in Sj4j In this restricted polynomial framework, the LP approach has interesting ad- 
ditional features that can be exploited for effective numerical computation. Indeed, 
what makes this LP approach attractive is that for the class of OCPs considered: 

• Only the moments of the o.m.'s appear in the LP formulation, so that we 
already end up with countably many variables, a significant progress. 

• Constraints on the support of the o.m.'s translate easily into either LP or SDP 
(Semi Definite Programming) necessary constraints on their moments. Even more, 
for (semi-algebraic) compact supports, relatively recent powerful results from real 
algebraic geometry make these constraints also sufficient. 

• When truncating to finitely many moments, the resulting LP or SDP's are 
solvable and their optimal values form a monotone nondecreasing sequence (indexed 
by the number of moments considered) of lower bounds on the optimal value of the 
LP (and thus of the OCP). 

Therefore, based on the above observations, we propose an approximation of the 
optimal value of the OCP via solving a hierarchy of SDPs (or linear matrix inequal- 
ities, LMIs)that provides a monotone nondecreasing sequence of lower bounds on 
the optimal value of the weak LP formulation of the OCP. In adddition, under some 
compactness assumption of the state and control constraint sets, the sequence of 
lower bounds is shown to converge to the optimal value of the LP, and thus the 
optimal value of the OCP when the former and latter are equal. 
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As such, it could be seen as a complement to the above shooting or direct meth- 
ods, and when the sequence of lower bounds converges to the optimal value, a test 
of their efficiency. Finally this approach can also be used to provide a certificate of 
unfeasibility. Indeed, if in the hierarchy of LMI-relaxations of the minimum time 
OCP, one is infeasible then the OCP itself is infeasible. It turns out that some- 
times this certificate is provided at an early stage in the hierarchy, i.e. with very 
few moments. This is illustrated on two simple examples. 

In a pioneering paper, Dawson jlO] had suggested the use of moments in the LP 
approach with o.m.'s, but results on the K-moment problem by Schmiidgen [38] and 
Putinar [37] were either not available at that time. Later, Helmes and Stockbridge 
[15] and Helmes, Rohl and Stockbridge [16] have used LP moment conditions for 
computing some exit time moments in some diffusion model, whereas for the same 
models, Lasserre and Prieto-Rumeau [29 have shown that SDP moment conditions 
are superior in terms of precision and number of moments to consider; in [30] . 
they have extended the moment approach for options pricing problems in some 
mathematical finance models. More recently, Lasserre, Prieur and Henrion |31j 
have used the o.m. approach for minimum time OCP without state constraint. 
Preliminary experimental results on Brockett's integrator example, and the double 
integrator show fast convergence with few moments. 

2. Occupation measures and the LP approach 

2.1. Definition of the optimal control problem. Let n and m be nonzero 
integers. Consider on K™ the control system 

(2.1) x(t) = f(t,x(t),u(t)), 

where / : [0, +oo) x R™ x K™ 1 — > K™ is smooth, and where the controls are bounded 
measurable functions, defined on intervals [0,T(u)] of K + , and taking their values 
in a compact subset U of W n . Let xq G W 1 , and let X and K be compact subsets 
of R™. For T > 0, a control u is said admissible on [0, T] whenever the solution x(-) 
of (|2.ip . such that x(0) = xq, is well defined on [0,T], and satisfies 

(2.2) (x(t),u(t)) e X x U a.e. on [0, T], 
and 

(2.3) x{T) e K. 

Denote by Ut the set of admissible controls on [0, T]. 

For u e Ut, the cost of the associated trajectory x(-) is defined by 

(2.4) J(0,T,x ,u) = [ h(t,x{t),u(t))dt + H(x(T)), 

Jo 

where h : [0, +oo) x R" x E m — ► M and H : W 1 -»■ M are smooth functions. 

Consider the optimal control problem of determining a trajectory solution of 
(|2.1I starting from x(0) = xq, satisfying the state and control constraints (|2.2|) . the 
terminal constraint (|2.3|) . and minimizing the cost (|2.4|) . The final time T may be 
fixed or not. 

If the final time T is fixed, we set 



(2.5) 



J*(0,7> ) := inf J{0,T,x o ,u), 
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and if T is free, we set 

(2.6) J*(0,x o ) := inf J(0, T, x , u), 

Note that a particular OCP is the minimal time problem from xo to K, by letting 
h = 1, H = 0. In this particular case, the minimal time is the first hitting time of 
the set K. 

It is possible to associate a stochastic or deterministic OCP with an abstract 
infinite dimensional linear programming (LP) problem P together with its dual 
P* (see for instance Hernandez-lerma and Lasserre [IB] for discrete-time Markov 
control problems, and Vinter [44], Hernandez et al. [23] for deterministic optimal 
control problems, as well as many references therein). We next describe this LP 
approach in the present context. 

2.2. Notations and definitions. For a topological space X, let M.(X) be the 
Banach space of finite signed Borel measures on X, equipped with the norm of 
total variation, and denote by Ai(X) + its positive cone, that is, the space of finite 
measures on X. Let G(X) be the Banach space of bounded continuous functions on 
X, equipped with the sup-norm. Notice that when X is compact Hausdorff, then 
M{X) ~ C(X)*, i.e., M(X) is the topological dual of C{X). 

Let £ := [0, T] x X, S := £ x U, and let Ci(£) be the Banach space of functions 
if G C(£) that are continuously differentiable. For ease of exposition we use the 
same notation g (resp. h) for a polynomial g G R[t, x] (resp. h G R[x]) and its 
restriction to the compact set £ (resp. K). 

Next, with u G U, let A : Ci(£) -> C(S) be the mapping 

(2.7) <p ^ A<p(t,x,u) := ^(t,x) + (f(t, x, u), V x¥ >(t, x)). 

Notice that d<p/dt+ {S7 x ip, /) G C(S) for every <p G Ci(S), because both X and U 
are compact, and / is understood as its restriction to S. 

Next, let C : Ci(E) -» C(S) x C(K) be the linear mapping 

(2.8) <p » Op := (-A<p,<pr), 

where <^t(x) := tp(T,x), for all x G X. Obviously, C is continuous with respect to 
the strong topologies of Ci(E) and C(S) x C(K). 

Both (C(S), jM(S)) and (C(K),7W(K)) form a rfuo( pozr of vector spaces, with 
duality brackets 

(h,/i) = Jhdn, V(/i,m) eC(S) xM(S), 

and 

( 5 , v)= J gdv, V ( 3 , G C(K) x M(K) 

Let £* : M(S) x M (K) d(E)* be the adjoint mapping of £, defined by 

(2.9) (M,&P) = (C*(n,v),<p), 
for all ((fx,v),ip) G M(S) x M(K) x C X (E). 

Remark 2.1. (i) The mapping is continuous with respect to the weak topolo- 
gies a(M(S) x M(K),C(S) x C(K)), and cr(Ci(E)*, C?i(E)). 
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(ii) Since the mapping C is continuous in the strong sense, it is also continu- 
ous with respect to the weak topologies cr(Ci(S), Ci(£)*) and <j(C(S) x 
C(K),M{S) x M(K)). 
(hi) In the case of a free terminal time T < To, it suffices to redefine C : 
Ci(S) -> C(S) x C([0,T ] x K) by £y> := (-A<p,<p). The operator £* : 
M(S)xM([0,T ]xK) -> Ci(S)* is still defined by f2T5|l . for all ((n,v),<p) £ 
M(S) x M([0,T ] x K) x Ci(S). 

For time-homogeneous free terminal time problems, one only needs func- 
tions <p of a; only, and so £ = S = X x U and C : Ci(£) -> C(S) x C(K). 

2.3. Occupation measures and primal LP formulation. Let T > 0, and let 

u = {u(t), < £ < T} be a control such that the solution of (|2.ip . with x(0) = xq, 
is well defined on [0, T]. Define the probability measure v u on R™, and the measure 
H u on [0,T] x E™ x M m , by 

(2.10) ^ U (.D) := I B [x(T)], DeB n , 

(2.11) (i u (ixBxC) := / luxe [(*(*), u(*))]*i 

J[o,T]nA 

for all rectangles (Ax B x C), with (A, £?, C) £ .4 x B„ x B m , and where £> n (resp. 
B m ) denotes the usual Borel er-algebra associated with 1" (resp. R m ), and A the 
Borel er-algebra of [0, T], and Ib(«) the indicator function of the set B. 

The measure fi u is called the occupation measure of the state-action (determin- 
istic) process (t,x(t),u(t)) up to time T, whereas v 11 is the occupation measure of 
the state x(T) at time T. 

Remark 2.2. If the control u is admissible on [0,T], i.e., if the trajectory x(-) 
satisfies the constraints (|2.2p and (|2.3|) , then v u is a probability measure supported 
on K (i.e. v u £ M(K)+), and fi u is supported on [0, T] x X x U (i.e. fi u £ -M(S)+). 
In particular, T = /Lt u (S). 

Conversely, if the support of /i u is contained in S = [0, T] x X x U and if 
^i u (S) = T, then (x(t),u(t)) £ X x U for almost every £ £ [0,T]. Indeed, with 

T = Ix.xu[(x(s),u(s))]ds 
Jo 

=>Ixxu[(z(s), «(«))] = 1 a.e. in [0,T], 
and hence (x(£), it(t)) £ X x U, for almost every £ £ [0, T]. If moreover the support 
of v u is contained in K, then x(T) £ K. Therefore, u is an admissible control on 
[0,T]. 

Then, observe that the optimization criterion (|2.5| now writes 

J(0,T> o ,u) = / Hdv u + [ hd[i u = ((/i u ,i/ u ), (&,#)}, 
Jk Js 

and one infers from (f2~Tj) , f[2T2]) and |2~3)l that 

(2-12) J^g T dv u -g(0,xo) = J s (|f + <V*$, />) 

for every <? £ Ci(S) (where gr(x) = g(T,x) for every a; £ K), or equivalently, in 
view of (EU) and p^]) . 

( 5 ,r(/i u ,i, u )} = v<?£d(s). 
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&p < (h,H) 



This in turn implies that 

Therefore, consider the infinite-dimensional linear program P 

(2.13) P: inf {((»,v),(h,H)) \ £*([/■, v) = ^(o,x )} 

(where A := Ai(S) + x A4(K) + ). Denote by inf P its optimal value and minP is 
the infimum is attained, in which case P is said to be solvable. The problem P is 
said feasible if there exists (/i, is) € A such that C*(p, v) = Sm !Xo ). 

Note that P is feasible whenever there exists an admissible control. 

The linear program P is a rephrasing of the OCP (|2.ip ~ (|2.5p in terms of the 
occupation measures of its trajectories (t, x(t),u(t)). Its dual LP reads 

(2.14) P*: sup {(S M ,ip)\ &p<(h,H)} 
where 

Aip(t, x, u) + h(t, x, u) >0 y(t, x,u)gS 

<p(T,x) <H(x) VieK 

Denote by supP* its optimal value and maxP* is the supremum is attained, i.e. if 
P* is solvable. 

Discrete-time stochastic analogues of the linear programs P and P* are also 
described in Hernandez-Lerma and Lasserre 18, 19J for discrete time Markov control 
problems. Similarly see Cho and Stockbridge [8], Kurtz and Stockbridge [27], and 
Hclmes and Stcokbridge [IS] for some continuous-time stochastic problems. 

Theorem 2.3. IfP is feasible, then: 

(i) P is solvable, i.e., inf P = minP < J(0, T, Xo). 

(ii) There is no duality gap, i.e., supP* = minP. 

(hi) // moreover, for every (t, x) G S, the set f(t, x, U) C R" is convex, and the 
function 

v >-> gt,x(v) ■= inf {h{t,x,u) : v = f(t,x,u)} 

is convex, then the OCP i2.1]) - [K5\) has an optimal solution and 
supP* = infP = minP = J*(0,T,x o ). 
For a proof see £15.41 Theorem I2.3f iii) is due to Vinter [33] . 

3. SEMIDEFINITE PROGRAMMING RELAXATIONS OF P 

The linear program P is infinite dimensional, and thus, not tractable as it 
stands. Therefore, we next present a relaxation scheme that provides a sequence of 
scmidcfinite programming, or linear matrix inequality relaxations (in short, LMI- 
relaxations) {Q r }, each with finitely many constraints and variables. 

Let R[x] = [xi, . . . x„] (resp. R[i, x, u] = R[t, xi, ■ . ■ x n , u%,..., u m \) denote the 
set of polynomial functions of the variable x (resp., of the variables t,x,u). 

Assume that X and K (resp., U) are compact semi- algebraic subsets of K" (resp. 
of M m ), of the form 



(3.1) X 

(3.2) K 

(3.3) U 



= {i£l" | v j (x)>0, j e J}, 
= {xeW n I 9 j (x)>0, jgJt}, 
= {ueR m | ^(u)>0, j e W}, 
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for some finite index sets Jt, J and W, where Vj, 9j and Wj are polynomial func- 
tions. Define 

(3.4) d(X,K,U):= max (deg# 7 -, deg vi, deg Wk)- 

jeJi, leJ, kew 

To highlight the main ideas, in this section we assume that /, h and H are 
polynomial functions, that is, h e R[t, x,u], H £ R[x], and / : [0, +oo) xR n xR m — > 
R" is polynomial, i.e., every component of / satisfies £ R[i, x, u], for fc = 1, . . . , n. 

3.1. The underlying idea. Observe the following important facts. 
The restriction of R[t, x] to £ belongs to Ci(S). Therefore, 

£*(/■»> v ) = <Wo) ^ (.9,^*(M^)) = 3(0,X ), Vj£f[t,4 

because S being compact, polynomial functions are dense in Ci(X) for the sup- 
norm. Indeed, on a compact set, one may simultaneously approximate a function 
and its (continuous) partial derivatives by a polynomial and its derivatives, uni- 
formly (see [24] pp. 65-66). Hence, the linear program P can be written 

P : { (/*,")£ A 

I s.t. (g,£*(fx,v)) = g(Q t x ), Vj6R[t,4 

or, equivalently, by linearity, 

3.5 P : ^ (p^)eA 

[ s.t. (£g,((i,v)} = g(Q,xo), Vg = (P x a ); (p, a) e N x N". 

The constraints of P, 

(3.6) <£ fls (/* ) i/)>= fl (0,zo) ) V F (f/);(p, a )eNxN», 

define countably many linear equality constraints linking the moments of /x and z/, 
because if g is polynomial then so are dg/dt and dg/dxk, for every and (V x g, /)■ 
And so, £3 is polynomial. 

The functions h,H being also polynomial, the cost ((/x, v), (h,H)) of the OCP 
(|2.ip ~ (|2.5[) is also a linear combination of the moments of /i and v. 

Therefore, the linear program P in (|3.5p can be formulated as a LP with count- 
ably many variables (the moments of fi and v) , and countably many linear equality 
constraints. However, it remains to express the fact that the variables should be 
moments of some measures fi and v, with support contained in S and K respectively. 

At this stage, one will make some (weak) additional assumptions on the poly- 
nomials that define the compact semi-algebraic sets X, K and U. Under such 
assumptions, one may then invoke recent results of real algebraic geometry on the 
representation of polynomials positive on a compact set, and get necessary and 
sufficient conditions on the variables of P to be indeed moments of two measures 
fi and v, with appropriate support. We will use Putinar's Positivstellensatz [37] 
described in the next section, which yields SDP constraints on the variables. 

One might also use other representation results like e.g. Krivine [55] and Vasilescu 
[43] , and obtain linear constraints on the variables (as opposed to SDP constraints) . 
This is the approach taken in e.g. Helmes et al. [16] . However, a comparison of the 
use of LP-constraints versus SDP constraints on a related problem [29] has dictated 
our choice of the former. 
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Finally, if g in (|3.6j) runs only over all monomials of degree less than r, one then 
obtains a corresponding relaxation Q r of P, which is now a finite-dimensional SDP 
that one may solve with public software packages. At last, one lets r — > oo. 

3.2. Notations, definitions and auxiliary results. For a multi-index a — 
(ai, . . . , a n ) £ N n , and for x = {x\, . . . , x n ) £ R n , denote x a :— x" 1 ■ ■ ■ x" n . Con- 
sider the canonical basis {a; Q } Qe N" (resp., {i p x Q M /3 } pe N,aeN",/3eN m ) of M.[x] (resp., 
of R[t,x,u]). 

Given two sequences y = {y a }aeN™ and z = {z 1 } 7£ nxN"xi™ of real numbers, 
define the linear functional L y : M[x] — > R by 

H(:= H a x a ) i ► L y {H) := ]T H a y a , 

ctGN™ qGN" 

and similarly, define the linear functional L z : M[i,a;,M] — > K by 

h^L z (h) := ^ h 1 z 1 = ^ h paf 3 z pa( 3, 

7GNxN™xN m pSN,aGN™,/3GN m 

where h(t,x,u) = E pG N,aeN" : /3eN™ h pa p t p x a u 13 . 

Note that, for a given measure v (resp., /i) on K (resp., on M x R ra x M m ), there 
holds, for every H £ M.[x] (resp., for every h £ M.[t, x, u]), 

(v, H)= [ Hdv= / V H a x a dv = V H aVa = L y (H), 

JR JR 

where the real numbers y a = J x a dv are the moments of the measure v (resp., 
(n, h) = L z (h), where z is the sequence of moments of the measure [i). 

Definition 3.1. For a given sequence z = {z T } 7 gNxN»xN™ of real numbers, the 
moment matrix M r {z) of order r associated with z, has its rows and columns 
indexed in the canonical basis {t p x a u@}, and is defined by 

(3.7) M P (z)( 7 ,/3) = z 7+/3 , x^Nxf x N m , \j\,\(3\<r, 

where I7I := jj. The moment matrix M r (y) of order r associated with a given 
sequence y = {y Q } Q gN™, has its rows and columns indexed in the canonical basis 
{x a }, and is defined in a similar fashion. 

Note that, if z has a representing measure //, i.e., if z is the sequence of moments 
of the measure \i on M x K" x R m , then L z {h) = J hdfj,, for every h £ R[t, x, u], and 
if h denotes the vector of coefficients of h £ R[t, x, u] of degree less than r, then 

(h,M r (z)h) = L z (h 2 ) = Jh 2 dn>0. 

This implies that M r (z) is symmetric nonnegative (denoted M r (z) >z 0), for every 
r. The same holds for M r (y). 

Conversely, not every sequence y that satisfies M r {y) y for every r, has a 
representing measure. However, several sufficient conditions exist, and in particular 
the following one, due to Berg [3] . 

Proposition 3.2. If y = {y a }aeN n satisfies \y a \ < 1 for every a £ N n , and 
M r (y) >: for every integer r, then y has a representing measure on W 1 , with 
support contained in the unit ball [—1,1]™. 

We next present another sufficient condition which is crucial in the proof of our 
main result. 
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Definition 3.3. For a given polynomial 8 <E R[f, x, it], written 

8(t,x,u)= O s t p x a u ff , 

define the localizing matrix M r (8 z) associated with z, 0, and with rows and columns 
also indexed in the canonical basis of R[t, x, u], by 

(3.8) M r (0*)(7,j9) = X)fc*a+7+/3 7,/JeNxrxN" 1 , \~/\,\P\<r. 

s 

The localizing matrix M r (8y) associated with a given sequence y = {y a }aeN n is 
defined similarly. 

Note that, if z has a representing measure /i on E x R™ x R m with support 
contained in the level set {(t, x, u) : 6(t, x, u) > 0}, and if h G R[t, x, u] has degree 
less than r, then 

(h,M r (6,z)h) = L z {8h 2 ) = J 8h 2 d^>0. 

Hence, M r (8 z) > 0, for every r. 

Let S 2 C M[x] be the convex cone generated in R[x] by all squares of polynomials, 
and let O C M. n be the compact basic semi-algebraic set defined by 

(3.9) Q:={xeR n \ gj(x)>0, j = l,...,m} 
for some family {gj}™= 1 C M.[x]. 

Definition 3.4. The set fi C 1" defined by (|3 ,9[) satisfies Putinar's condition if 
there exists u S M.[x] such that u = uq + u j9j f° r some family {ujY" =Q C X 2 , 

and the level set {x £ R™ | u(x) > 0} is compact. 

Putinar's condition is satisfied if e.g. the level set {x : gu{x) > 0} is compact 
for some k, or if all the gj's are linear, in which case SI is a polytope. In addition, 
if one knows some M such that < M whenever x € SI, then it suffices to add 
the redundant quadratic constraint M 2 — \\x\\ 2 > in the definition Q3.9P of SI, and 
Putinar's condition is satisfied (take u := M 2 — ||a;|| 2 ). 

Theorem 3.5 (Putinar's Positivstellensatz [37]). Assume that the set fl defined by 
(3. Sty satisfies Putinar's condition. 

(a) If f S E[x] and f > on SI, then 

m 

(3-10) / = /o + 

3=1 

for some family {fj}™= C S 2 . 

(b) Let y = {y a }a£N n be a sequence of real numbers. If 

(3.11) M r (y)t0; M r ( 9j y) h0, j = l,...,m; Vr = 0,l,... 

then y has a representing measure with support contained in SI. 
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3.3. LMI-relaxations of P. Consider the linear program P defined by (|3.5|) . 

Since the semi- algebraic sets X, K and U defined respectively by (|3.1|) . ()3.2j) 
and (|3.3p are compact, with no loss of generality, we assume (up to a scaling of the 
variables x, u and i) that T = 1, X,K C [-1, 1]™ and U C [-1, l] m . 



Next, given a sequence z = {z 7 } indexed in the basis of 



denote z(t), 



z(x) and z(u) its marginals with respect to the variables t, x and u, respectively. 
These sequences are indexed in the canonical basis of R[t] , M[x\ and M.[u] repectively. 
For instance, writing 7 = (k, a, (3) G N x N" x N n , 



{ z (t)} = {^fc,0,o}fcGN; {z{x)} = {z ,Q,o}ae 



{z(u)} = {z o . O:0 } 0i 



(3.12) Q r 



Let ro be an integer such that 2ro > max (deg /, deg h, deg H, 2d(X, K,U)), 
where d(X, K, U) is defined by (|3.4|) . For every r > ro, consider the LMI-relaxation 

inf LJh)+LJH) 
M r (y), M r (z) y 

M r _do g fi 3 -(0jj/) b 0, jeJi 

M r _de g! ,3 (fj «(a;)) h 0, je J 

Mr-dcg«, fc (Wfe z(t<)) b 0, k £ W 

M r _i(t(l-i)-r(t)) >r 

L y ( 9l ) - L,(9 5 /3t + (V x <7, /}) - ,g(0, aro), V 5 = 
with p + |a| — 1 + deg / < 2r 

whose optimal value is denoted by inf Q r . 

OCP with free terminal time. For the OCP (|2.6p . i.e., with free terminal time 
T < To, we need adapt the notation because T is now a variable. As already men- 
tioned in Remark l2.1f iii). the measure v in the infinite dimensional linear program 
P defined in (|2. 13|) . is now supported in [0, To] x K (and [0,1] x K after re-scaling) 
instead of K previously. Hence, the sequence y associated with v is now indexed in 
the basis {t p x a } of M[t, x] instead of {x a } previously. Therefore, given y — {yka} 
indexed in that basis, let y(t) and y(x) be the subsequences of y defined by: 

y(t) := {y k o}k, k G N; ; y{x) = {y 0c[ }, a G N n . 

Then again (after rescaling), the LMI-relaxation Q r now reads 

L z (h) + L y {H) 

M r (z) y 



(3.13) 



inf 

y,z 

M r [y) 

M, 
M, 



j G Ji 

r(vj){Vj z(x)) h 0, jeJ 

r(w k ){Wk z{u)) >: 0, k G W 



M^CL-i)^)) h0 
M r _x(t(l -t)z(t)) h 
L y {g)-L z (dg/dt+(V x gJ)) 
with p + \a\ — 1 + deg / < 2r 



g(0,x o ), Vg = (tPx a ) 



The particular case of minimal time problem is obtained with h = 1, H = 0. 

For time-homogeneous problems, i.e., when h and / do not depend on t, one may 
take (resp. v) supported on X x U (resp. K), which simplifies the associated 
LMI-relaxation ([37TE]) . 

The main result is the following. 
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Theorem 3.6. Let X, K C [—1,1]", and U C [— l,l] m be compact basic semi- 
algebraic-sets respectively defined by 13. VS. 2i) and \3. 3\) . Assume that X, K 
and U satisfy Putinar's condition (see Definition \3.4\l ), and let Q r be the LMI- 
relaxation defined in H3.12\) . Then, 

(i) inf Q r f minP as r —> oo; 

(ii) if moreover, for every (t,x) G S, the set f(t,x,XJ) C K™ is convex, and the 
function 

v h-> fjt.x(v) := inf { h(t, x,u) \ v = f(t, x, u)} 

is convex, then inf Q r f minP = J* (0,T, Xq), as r — > oo. 
The proof of this result is postponed to the Appendix in Section §5.51 

3.4. The dual Q*. We describe the dual of the LMI-relaxation Q r which is also 
a semidefinite program, denoted Q*, and relate Q* with the dual P* of P, defined 
in (|2T4| . 

Let s(r) be the cardinal of the set V r := {(k, a) G N x N n | k + \a\ < r - r }, 
and given A G R s(r) , let A r G be the polynomial 



(i, a;) i— > A r (i, x) := 



E 

(fe,a)ev r 



Afeo, t k x° 



Consider the semidefinite program: 

sup A r (0,x o ), 

h + A A r = 9o i(l - i) + Efcew 9fe ^ + Eje j 3j ty, 

,;-!. |.4) Q; .{ H- A r (l, .) = io + E,- 6 Ji «j ^ 

9o e R[t], ql e R[u], q] G R[x], l 3 G K[x] 

{qo, , Zo, s.o.s. (sums of squares polynomials), and 
degljOj, degqfVj, degq%w k , degq < 2r - 2; deg/ < 2r. 

The LMI Q* is a reinforcement of P* in the following sense: 

• the unknown function ip G Ci(S) is now replaced with a polynomial A r G 
M.[t, x] of degree less than 2r; 

• the constraint —Atp < /i for (t, x, u) G S, is now replaced with the constraint 
h + AA r > on S and the polynomial h + AA r > which is nonnegative 
on S, has Putinar's representation go t(l —t) + Efcew Ik Wk + Ejej Qj V J' 

• the constraint ip\ < H for x G K, is replaced with the constraint H — 
A r (l, •) > on K, and the polynomial H — A r (l, .) which is nonnegative 
on K, has Putinar's representation Iq + Eje./i lj Oj. 

Assume that Q* is solvable. A natural question is to know whether or not we 
can use an optimal solution 

QO, Qjilki fo) ij'i Ar 
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of Q* to obtain some information on an optimal solution of P. The most natural 
idea is to look for the zero set in S of the polynomial 



Indeed, under the assumptions of Theorem l3.6l if sup Q* = inf Q,., then A r (0, xq) f=a 
inf Q r « minP = supP*, and so, the polynomial A r £ R[i,x] seems to be a good 
candidate to approximate a nearly optimal solution (p £ Ci(S) ofP*. 

Next, as Q* is an approximation of a weak formulation of the HJB optimality 
equation, one may hope that the zero set in S of the polynomial h + AA r provides 
some good information on the possible states x*(t) and controls u*(t) at time t in 
an optimal solution of the OCP JUT]) -JUS]). 

That is, fixing an arbitrary t$ £ [0, 1], one may solve the equation 



and look for solutions (x, u) £ X x U. 

All these issues deserve further investigation beyond the scope of the present 
paper. However, at least in the minimum time problem for the (state and control 
constrained) double integrator example considered in §5.1( we already have some 
numerical support for the above claims. 

3.5. Certificates of non controllability. For minimum time OCPs, i.e., with 
free terminal time T and instantaneous cost h = 1, and H = 0, the LMI-relaxations 
Q r defined in (I3.13P may provide certificates of non controllability. 

Indeed, if for a given initial state xq £ X, some LMI-relaxation Q r in the hierar- 
chy has no feasible solution, then the initial state xo cannot be steered to the origin 
in finite time. In other words, inf Q r = +oo provides a certificate of uncontrolla- 
bility of the initial state xo. It turns out that sometimes such certificates can be 
provided at cheap cost, i.e., with LMI-relaxations of low order r. This is illustrated 
on the Zermelo problem in £ 15.31 

Moreover, one may also consider controllability in given finite time T, that is, 
consider the LMI-relaxations as defined in (|3 . 1 2[) with T fixed, and H = 0, h = 1. 
Again, if for a given initial state xo £ X, the LMI-relaxation Q r has no feasible 
solution, the initial state xo cannot be steered to the origin in less than T units of 
time. And so, inf Q r = +oo also provides a certificate of uncontrollability of the 
initial state Xo- 

4. Generalization to smooth optimal control problems 

In the previous section, we assumed, to highlight the main ideas, that /, h and 
H were polynomials. In this section, we generalize Theorem l3.6l and simply assume 
that /, h and H are smooth. Consider the linear program P defined in the previous 
section 



Since the sets X, K and U, defined previously, are compact, it follows from [9] 
(see also [2H PP- 65-66]) that / (resp. h, resp. H) is the limit in Ci(S) (resp. Ci(S), 
resp. Ci(K)) of a sequence of polynomials f p (resp. h p , resp. H p ) of degree p, as 
p — > +oo. 



(t, x, u) i * q t(l -t)+Y^ 1k w k + tf v i- 



kew je.J 



X! Qk( u ) w k(u) +^2qj(x)vj(x) = -q t (l-t ) 




s.t. (g,£*(n,v)) = g(0, x ), VjeR[t,i 
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Hence, for every integer p, consider the linear program P p 

f inf {{{^v),{h p ,H p )) 
p : ) (M^)eA 

[ s.t. (g,C;(jt,v)) = g(0,x o ), Vj6R[i,i], 

the smooth analogue of P, where the linear mapping C p : Ci(S) — > C(S) x C(K) 
is defined by 

Cpf := (-A P ¥J,^ T ), 
and where A p : Ci(S) — > C(S) is defined by 

A p <,5(t, a;, u) := -^{t^) + (fp(t,x,u), W x ip{t,x)). 

For every integer r > max(p/2, d(X, K, U)), let Q r , p denote the LMI-relaxation 
(|3.12p associated with the linear program P p . 

Recall that from Theorem 13.61 if K, X and U satisfy Putinar's condition, then 
inf Qr )P t minPp as r — > +00; 

The next result, generalizing Theorem l3.61 shows that it is possible to let p tend 
to +00. For convenience, set 

Vr lP = inf Q r ,p, v p = minPp, v = minP. 

Theorem 4.1. Let X, K c [— 1, 1]", and U C [—1, l] m be compact semi- algebraic- 
sets respectively defined by \3.1\) , US. 2\) and VS. 3\) . Assume that X, K and U satisfy 
Putinar's condition (see Definition \3.J$ ). Then, 

(i) v= lim lim v r p = lim sup v r p < J*(0,T,Xq). 

p— >+oo r— >+oo ' p-i+oo . /9 
2r>p r>p/z 

(ii) Moreover if for every (t,x) G S, the set f(t,x,XJ) C K™ is convex, and the 
function 

v i-> gt,x(v) := inf { h(t, x,u) \v = f(t, x, u)} 

is convex, then v — J*(0,T,xq)- 
The proof of this result is in the Appendix, Section i j5.6l 

From the numerical point of view, depending on the functions /, h, H, the 
degree of the polynomials of the approximate OCP P p may be required to be large, 
and hence the hierarchy of LMI-relaxations (Q r ) in (|3 . 1 2[) might not be efficiently 
implcmentable, at least in view of the performances of public SDP solvers available 
at present. 

Remark 4.2. The previous construction extends to smooth optimal control problems 
on Riemannian manifolds, as follows. Let M and N be smooth Riemannian mani- 
folds. Consider on M the control system (|2~Tj) . where / : [0, +00) xMxN — ► TM 
is smooth, and where the controls are bounded measurable functions, defined on 
intervals [0, T(u)] of M + , and taking their values in a compact subset U of N. Let 
x S M, and let X and K be compact subsets of M. Admissible controls are de- 
fined as previously. For an admissible control u on [0,T], the cost of the associated 
trajectory x(-) is defined by (|2.4p . where where h : [0, +00) x M x N — > K and 
H : M — > K are smooth functions. 

According to Nash embedding Theorem t 33] , there exist an integer n (resp. m) 
such that M (resp. N) is smoothly isometrically embedded in R™ (resp. R m ). In 
this context, all previous results apply. 
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This remark is important for the applicability of the method described in this 
article. Indeed, many practical control problems (in particular, in mechanics) are 
expressed on manifolds, and since the optimal control problem investigated here is 
global, they cannot be expressed in general as control systems in M. n (in a global 
chart) . 

5. Illustrative examples 

We consider here the minimal time OCP, that is, we aim to approximate the 
minimal time to steer a given initial condition to the origin. We have tested the 
above methodology on two test OCPs, the double and Brockett integrators, for 
which the associated optimal value T* can be calculated exactly. The numerical 
examples in this section were processed with our Matlab package GloptiPoly 50. 

5.1. The double integrator. Consider the double integrator system in R 2 
1 j x 2 (t) = u(t), 

where x — (x\, x 2 ) is the state and the control u = u(t) £ U, satisfies the constraint 
\u{t)\ < 1, for all t > 0. In addition, the state is constrained to satisfy x 2 {t) > —1, 
for all t. In this time-homogeneous case, and with the notation of Section [21 we 
have X = {x e M 2 : x 2 > -1}, K = {(0, 0)}, and U = [-1, 1]. 

Remark 5.1. The theorem obviously extends, up to scaling, to the case of arbitrary 
compact subsets X, K C R™ and U C [-1, l] m . 

Observe that X is not compact and so the convergence result of Theorem 13.61 
may not hold. In fact, we may impose the additional constraint Hx^Hoo < M for 
some large M (and modify X accordingly), because for initial states xq with ||xo||oo 
relatively small with respect to M, the optimal trajectory remains in X. How- 
ever, in the numerical experiments, we have not enforced an additional constraint. 
We have maintained the original constraint x 2 > — 1 in the localizing constraint 
M r _ r ( v .)(vjz(x)) y 0, with x i — ► Vj(x) = x 2 + 1. 

5.1.1. Exact computation. For this very simple system, one is able to compute 
exactly the optimal minimum time. Denoting T{x) the minimal time to reach the 
origin from x = {x\, x 2 ), we have: 

If x-i > l-xJ/2 and x 2 > -1 then T(x) = xI/2+x-l+x 2 + 1. If -^/2signx 2 < 
Xi < 1 — x^/2 and x 2 > —1 then T(x) — 2\Jx\j2 + x\ + x 2 . If x\ < —x\j2 signa^ 
and x 2 > — 1 then T{x) = 2y 'x\/2 — x\ — x 2 . 

5.1.2. Numerical approximation. Table [1] displays the values of the initial state 
x Q G X, and denoting inf Q r (xo) the optimal value of the LMI-relaxation (|3.13[) for 
the minimum time OCP (|5.ip with initial state xq, Tables [U [31 and [H display the 
numerical values of the ratii inf Q r (xo)/T(xo) for r = 2, 3 and 5 respectively. 

In Figures [U [21 and [3] one displays the level sets of the ratii inf Q r /T(xo) for 
r = 2, 3 and 5 respectively. The closer to white the color, the closer to 1 the ratio 
inf Q r /T(xo). 

Finally, for this double integrator example we have plotted in Figure 2] the level 
sets of the function Ag(x) — T(x) where T(x) is the known optimal minimum time 



GloptiPoly 3 can be downloaded at www.laas.fr/~henrion/software 
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0.0 


0.2 


0.4 


0.6 


0.8 


1.0 


1.2 


1.4 


1.6 


1.8 


2.0 


X02 


-1.0 


-0.8 


-0.6 


-0.4 


-0.2 


0.0 


0.2 


0.4 


0.6 


0.8 


1.0 



Table 2. Double integrator: ratio inf Q2/T(x ) 



second LMI-rclaxation: r=2 



0.4598 


0.3964 


0.3512 


0.9817 


0.9674 


0.9634 


0.9628 


0.9608 


0.9600 


0.9596 


0.9595 


0.4534 


0.3679 


0.9653 


0.9347 


0.9355 


0.9383 


0.9385 


0.9386 


0.9413 


0.9432 


0.9445 


0.4390 


0.9722 


0.8653 


0.8457 


0.8518 


0.8639 


0.8720 


0.8848 


0.8862 


0.8983 


0.9015 


0.4301 


0.7698 


0.7057 


0.7050 


0.7286 


0.7542 


0.7752 


0.7964 


0.8085 


0.8187 


0.8351 


0.4212 


0.4919 


0.5039 


0.5422 


0.5833 


0.6230 


0.6613 


0.6870 


0.7121 


0.7329 


0.7513 


0.0000 


0.2238 


0.3165 


0.3877 


0.4476 


0.5005 


0.5460 


0.5839 


0.6158 


0.6434 


0.6671 


0.4501 


0.3536 


0.3950 


0.4403 


0.4846 


0.5276 


0.5663 


0.5934 


0.6204 


0.6474 


0.6667 


0.4878 


0.4493 


0.4699 


0.5025 


0.5342 


0.5691 


0.5981 


0.6219 


0.6446 


0.6647 


0.6824 


0.5248 


0.5142 


0.5355 


0.5591 


0.5840 


0.6124 


0.6312 


0.6544 


0.6689 


0.6869 


0.7005 


0.5629 


0.5673 


0.5842 


0.6044 


0.6296 


0.6465 


0.6674 


0.6829 


0.6906 


0.7083 


0.7204 


0.6001 


0.6099 


0.6245 


0.6470 


0.6617 


0.6792 


0.6972 


0.7028 


0.7153 


0.7279 


0.7369 



Table 3. Double integrator: ratio inf Q3/T(xq) 



third LMI-rclaxation: r=3 


0.5418 


0.4400 


0.3630 


0.9989 


0.9987 


0.9987 


0.9985 


0.9984 


0.9983 


0.9984 


0.9984 


0.5115 


0.3864 


0.9803 


0.9648 


0.9687 


0.9726 


0.9756 


0.9778 


0.9798 


0.9815 


0.9829 


0.4848 


0.9793 


0.8877 


0.8745 


0.8847 


0.8997 


0.9110 


0.9208 


0.9277 


0.9339 


0.9385 


0.4613 


0.7899 


0.7321 


0.7401 


0.7636 


0.7915 


0.8147 


0.8339 


0.8484 


0.8605 


0.8714 


0.4359 


0.5179 


0.5361 


0.5772 


0.6205 


0.6629 


0.7013 


0.7302 


0.7540 


0.7711 


0.7891 


0.0000 


0.2458 


0.3496 


0.4273 


0.4979 


0.5571 


0.5978 


0.6409 


0.6719 


0.6925 


0.7254 


0.4556 


0.3740 


0.4242 


0.4789 


0.5253 


0.5767 


0.6166 


0.6437 


0.6807 


0.6972 


0.7342 


0.4978 


0.4709 


0.5020 


0.5393 


0.5784 


0.6179 


0.6477 


0.6776 


0.6976 


0.7192 


0.7347 


0.5396 


0.5395 


0.5638 


0.5955 


0.6314 


0.6600 


0.6856 


0.7089 


0.7269 


0.7438 


0.7555 


0.5823 


0.5946 


0.6190 


0.6453 


0.6703 


0.7019 


0.7177 


0.7382 


0.7539 


0.7662 


0.7767 


0.6255 


0.6434 


0.6656 


0.6903 


0.7193 


0.7326 


0.7543 


0.7649 


0.7776 


0.7917 


0.8012 



Table 4. Double integrator: ratio inf Q5/T(:co) 



fifth LMI-rclaxation: r=5 


0.7550 


0.5539 


0.3928 


0.9995 


0.9995 


0.9995 


0.9994 


0.9992 


0.9988 


0.9985 


0.9984 


0.6799 


0.4354 


0.9828 


0.9794 


0.9896 


0.9923 


0.9917 


0.9919 


0.9923 


0.9923 


0.9938 


0.6062 


0.9805 


0.9314 


0.9462 


0.9706 


0.9836 


0.9853 


0.9847 


0.9848 


0.9862 


0.9871 


0.5368 


0.8422 


0.8550 


0.8911 


0.9394 


0.9599 


0.9684 


0.9741 


0.9727 


0.9793 


0.9776 


0.4713 


0.6417 


0.7334 


0.8186 


0.8622 


0.9154 


0.9448 


0.9501 


0.9505 


0.9665 


0.9637 


0.0000 


0.4184 


0.5962 


0.7144 


0.8053 


0.8825 


0.9044 


0.9210 


0.9320 


0.9544 


0.9534 


0.4742 


0.5068 


0.6224 


0.7239 


0.7988 


0.8726 


0.8860 


0.9097 


0.9263 


0.9475 


0.9580 


0.5410 


0.6003 


0.6988 


0.7585 


0.8236 


0.8860 


0.9128 


0.9257 


0.9358 


0.9452 


0.9528 


0.6106 


0.6826 


0.7416 


0.8125 


0.8725 


0.9241 


0.9305 


0.9375 


0.9507 


0.9567 


0.9604 


0.6864 


0.7330 


0.7979 


0.8588 


0.9183 


0.9473 


0.9481 


0.9480 


0.9559 


0.9634 


0.9733 


0.7462 


0.8032 


0.8564 


0.9138 


0.9394 


0.9610 


0.9678 


0.9678 


0.9696 


0.9755 


0.9764 




LMI 3 - Relative gep in % 




Figure 2. Double integrator: level sets inf Q 3 /T(a;o) 



to steer the initial state x to 0; the problem being time-homogeneous, one may take 
A r G M[x] instead of R[t, x] . For instance, one may verify that when the initial state 
is in the region where the approximation is good, then the whole optimal trajectory 
also lies in that region. 

5.2. The Brockett integrator. Consider the so-called Brockett system in K 3 

Xi(t) = Ui(t), 

(5.2) x 2 (t) = u 2 (t), 

X3{t) = u 1 (t)x 2 (t) - u 2 (t)x 1 (t), 
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Figure 3. Double integrator: level sets inf Q 5 /T(a;o) 




-2 -1.5 -1 -0.5 0.5 1 1.5 2 



Figure 4. Double integrator: level sets A$(x) — T{x) and optimal 
trajectory starting from xi(0) = ^(O) = 1. 

where x — {x\, x%, X3), and the control u = (ui(t),U2(t)) G hi, satisfies the con- 
straint 

(5.3) ui(t) 2 + u 2 (t) 2 < 1, Vt>0. 

In this case, we have X = M 3 , K = {(0, 0, 0)}, and U is the closed unit ball of R 2 , 
centered at the origin. 

Note that set X is not compact and so the convergence result of Theorem 13.61 
may not hold, see the discussion at the beginning of example 15. 11 Nevertheless, in 
the numerical examples, we have not enforced additional state constraints. 
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5.2.1. Exact computation. Let T(x) be the minimum time needed to steer an initial 
condition x € M 3 to the origin. We recall the following result of [2] (in fact given 
for equivalent (reachability) OCP of steering the origin to a given point x). 

Proposition 5.2. Consider the minimum time OCP for the system &5.ty) with 
control constraint \5. 3\) . The minimum time T(x) needed to steer the origin to a 
point x = {xi,X2,xz) S R 3 is given by 



^ 77 x 6^/xJ+xJ+2\xJ\ 

(5.4) T(xi,x 2 ,x 3 ) = 



\J 9 + sin 2 9 - sin 9 cos ( 



where 9 — 6(x±, x 2 , X3) is the unique solution in [0, 7r) of 

9 — sm9cos9, 7 2 . „. . 

(5.5) r^- (x 2 +x 2 ) = 2|x 3 |. 

sin 9 

Moreover, the function T is continuous on M 3 , and is analytic outside the line 
x\ = x 2 = 0. 

Remark 5.3. Along the line x\ = X2 — 0, one has 



T(0,0,x 3 ) = ^2ir\x 3 \. 

The singular set of the function T, i.e. the set where T is not C , is the line x\ = 
X2 = in ]R 3 . More precisely, the gradients dT/dxi, i = 1,2, are discontinuous at 
every point (0, 0, X3), x 3 ^ 0. For the interested reader, the level sets {(x±, x 2 , x 3 ) 6 
M 3 I T(xi,x 2 ,x 3 ) = r}, with r > 0, are displayed, e.g., in Prieur and Trelat [55] , 

5.2.2. Numerical approximation. Recall that the convergence result of Theorem 
13.61 is guaranteed for X compact only. However, in the present case X = K 3 is not 
compact. One possibility is to take for X a large ball of M 3 centered at the origin 
because for initial states xq with norm ||xo|| relatively small, the optimal trajectory 
remains in X. However, in the numerical experiments presented below, we have 
chosen to maintain X = R 3 , that is, the LMI-relaxation Q r does not include any 
localizing constraint M r _ r t v .\(vj z(x)) ^ on the moment sequence z(x). 

Recall that infQ r f minP as r increases, i.e., the more moments we consider, 
the closer to the exact value we get. 

In Table [5] we have displayed the optimal values inf Q r for 16 different values of 
the initial state x(0) = xq, in fact, all 16 combinations of xqi = 0, X02 = 0, 1, 2, 3, 
and X03 = 0, 1,2,3. So, the entry (2,3) of Table [5] for the second LMI-relaxation 
is inf Q2 for the initial condition xo = (0, 1, 2). At some (few) places in the table, 
the * indicates that the SDP solver encountered some numerical problems, which 
explains why one finds a lower bound inf Q r _i slightly higher than infQ r , when 
practically equal to the exact value T*. 

Notice that the upper triangular part (i.e., when both first coordinates Xo2,x 3 
of the initial condition are away from zero) displays very good approximations with 
low order moments. In addition, the further the coordinates from zero, the best. 

For another set of initial conditions xoi = 1 and xoi = {1,2,3} Table |6] displays 
the results obtained at the LMI-relaxation Q4. 

The regularity property of the minimal-time function seems to be an important 
topic of further investigation. 



OPTIMAL CONTROL AND LMI-RELAXATIONS 

Table 5. Brockett integrator: LMI-relaxations: inf Q r 



19 



first LMI-relaxation: r=l 


0.0000 


0.9999 


1.9999 


2.9999 


0.0140 


1.0017 


2.0010 


3.0006 


0.0243 


1.0032 


2.0017 


3.0024 


0.0295 


1.0101 


2.0034 


3.0040 


Second LMI-relaxation: r=2 


0.0000 


0.9998 


1.9997* 


2.9994* 


0.2012 


1.1199 


2.0762 


3.0453 


0.3738 


1.2003 


2.1631 


3.1304 


0.4946 


1.3467 


2.2417 


3.1943 


Third LMI-relaxation: r=3 


0.0000 


0.9995 


1.9987* 


2.9984* 


0.7665 


1.3350 


2.1563 


3.0530 


1.0826 


1.7574 


2.4172 


3.2036 


1.3804 


2.0398 


2.6797 


3.4077 


Fourth LMI-relaxation: r=4 


0.0000 


0.9992 


1.9977 


2.9952 


1.2554 


1.5925 


2.1699 


3.0478 


1.9962 


2.1871 


2.5601 


3.1977 


2.7006 


2.7390 


2.9894 


3.4254 


Optimal time T* 


0.0000 


1.0000 


2.0000 


3.0000 


2.5066 


1.7841 


2.1735 


3.0547 


3.5449 


2.6831 


2.5819 


3.2088 


4.3416 


3.4328 


3.0708 


3.4392 



Table 6. Brockett integrator: inf Q4 with xqi = 1 



fourth LMI-relaxation: r=4 


1.7979 


2.3614 


3.2004 


2.3691 


2.6780 


3.3341 


2.8875 


3.0654 


3.5337 


Optimal time T* 


1.8257 


2.3636 


3.2091 


2.5231 


2.6856 


3.3426 


3.1895 


3.1008 


3.5456 



5.3. Certificate of uncontrollabilty in finite time. Consider the so-called Zer- 
melo problem in K 2 studied in Bokanowski et al. [5] 



(5.6) 



x\{t) = 1 — a X2(t) + v cos 9 
X2(t) = v sin 9 
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FIGURE 5. Zermelo problem: uncontrollable states with Qi 

with a = 0.1. The state x is constrained to remain in the set X := [—6, 2] x [—2, 2] C 
M 2 , and we also have the control constraints < v < 0.44, as well as 9 £ [0, 2tt]. 
The target K to reach from an initial state xq is the ball B(0, p) with p :— 0.44. 

With the change of variable u\ ~ vcos9, U2 — vsuiQ, and U := {u : u\ + u\ < 
p 2 }, this problem is formulated as a minimum time OCP with associated hierarchy 
of LMI-relaxations (Q r ) defined in (|3 . 13[) . Therefore, if an LMI-relaxation Q r at 
some stage r of the hierarchy is infeasible then the OCP itself is infeasible, i.e., the 
initial state xq cannot be steered to the target K while the trajectory remains in 
X. 

Figures [5] and [6] display the uncontrollable initial states xq S X found with 
the infeasibility of the LMI-relaxations Qi and Q2 respectively. With Qi, i.e. 
by using only second moments, we already have a very good approximation of 
the controllable set displayed in [5], and Q2 brings only a small additional set of 
uncontrollable states. 
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X 2 



FIGURE 6. Zermelo problem: uncontrollable states with Q2 
Appendix 

5.4. Proof of Theorem 12.31 We first prove Item (i). Consider the linear program 
P defined in (|2.13|) . assumed to be feasible. From the constraint C*(ji, v) = <5(o,x )> 
one has 

J g{T, x)dv - J \ ^L(t, x) + (^(t, x)J(t, x, u))j = 5 (0, ar ), V 5 E 

In particular, taking g(t,x) = 1 and g(t,x) = T — t, it follows that fi(S) — T and 
f(K) = 1. Hence, for every (/i, v) satisfying C*{}i,v) = S^ Qxo ^, the pair (A/x, u) 
belongs to the unit ball B\ of (-M(S) x A4(K)). From Banach-Alaoglu Theorem, 
B\ is compact for the weak * topology, and even sequentially compact because B\ 
is metrizable (see e.g. Hernandez-Lerma and Lasserre [201 Lemma 1.3.2]). Since C* 
is continuous (see Remark l2.ip . it follows that the set of (//, v) satisfying £*(//, v) = 
%),x ) i s a dosed subset of Bi<~)(A4(S) + x A4(K) + ), and thus is compact. Moreover, 
since the linear program P is feasible, this set is nonempty. Finally, since the linear 
functional to be minimized is continuous, P is solvable. 
We next prove Item (ii). Consider the set 

D :={(C*(ii,u),((h,H),0i,u))) I (p,u) G M(S)+ x M(K) + }. 

To prove that D is closed, consider a sequence {(/x n , f n )}ti6ti of M(S) + x M(K)+ 
such that 

(5.7) {£*(» a ,v a ),((h,H),(fM a ,v a )}) -» (a, 6), 

for some (a, 6) G Ci(E)* xR. It means that £*(/i n , i/ n ) - * a, an d ((h, H), (fi n , v n )) — > 
6. In particular, taking tfQ := T — t and ipi = 1, there must hold 

//„(S) = (<po,£*(n n ,is n )) -> (930, a), ^n(K) = (ip\,C*(p, n ,u n )) —> (<pi,a). 

Hence, there exist no G N and a ball i? r of .M(S) x .M(K), such that (u n , z/„) G B r 
for every n > uq. Since £? r is compact, up to a subsequence (/i n ,^„) converges 
weakly to some (u, v) G .A4(S) + x A4(K) + . This fact, combined with (|5.7[) and the 
continuity of £*, yields a = £*(/x, ^), and b — ((h, H), (/1, u)). Therefore, the set D 
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is closed. 

From Anderson and Nash [TJ Theorems 3.10 and 3.22], it follows that there is 
no duality gap between P and P*, and hence, with (i), supP* = minP. 

Item (iii) follows from Vinter [44J Theorems 2.1 and 2.3], applied to the mappings 

F(t, x) := f(t, x, U) , l(t, x, v) :— inf { h(t, x,u) \v = /(t, x, u) }, 
for (t,x) G R x E". □ 

5.5. Proof of Theorem 13.61 First of all, observe that Q r has a feasible solution. 
Indeed, it suffices to consider the sequences y and z consisting of the moments (up 
to order 2r) of the occupations measures v u and /i u associated with an admissible 
control uel/of the OCP (j2~Tj) - (|2~5l) . 

Next, observe that, for every couple (y, z) satisfying all constraints of Q r , there 
must holds yo = 1 and Zq = 1. Indeed, it suffices to choose g(t,x) = 1 and 
g(t,x) = 1 — t (or t) in the constraint 

L y ( 9l ) - L z {dg/dt + (V x g, /)) = g{0, x ). 

We next prove that, for r sufficiently large, one has |z(a;) Q | < 1, |^(u)/3| < 1, 
|z(i)fc| < 1, for every k, and \y a \ < 1. We only provide the details of the proof for 
z(x), the arguments being similar for z(u), z(t) and y. 

Let E 2 C M.[x] be the space of sums of squares (s.o.s.) polynomials, and let 
Q C M.[x] be the quadratic modulus generated by the polynomials Vj £ R[x] that 
define X, i.e., 

Q := { cr e R[x] | a = (Ta + ^ajVj with aj € S 2 , V j e {0} U J}. 

In addition, let Q(t) C Q be the set of elements a of Q which have a representation 
o~o+^2j G j o~j Vj for some s.o.s. family {cj} C E 2 with degco £ 2i and degcrjUj < 2t 
for every j G J. 

Let r € N be fixed. Since X C [—1,1]", there holds 1 ± x a > on X, for every 
a G N n with \a\ < 2r. Therefore, since X satisfies Putinar' condition (see Definition 
I3.4p . the polynomial x i— > 1 ± x a belongs to Q (see Putinar [37]). Moreover, there 
exists l(r) such that x i — > 1 i x a G Q(l(r)) for every |a| < 2r. Of course, cc 
1 ± x a G <5(0 for every |a| < 2r, whenever I > i(r). 

For every feasible solution z of Q/( r ) one has 

\z{x) a \ = \L z {x a )\ <z = l, \a\ <2r. 
This follows from zq = 1, M/f r ^(z) >z and Mi^_ r ^ Vj ^(vj z(x)) >z 0, which implies 

m 

zq ± z(a;) a = L z (l ± x a ) = L z (a ) + }] L z ^ (<Tj Vj) > 0. 

3=1 

With similar arguments, one redefines l(r) so that, for every couple (y, z) satisfying 
the contraints of Q/( r ), one has 

< z k {t) < 1 and |a!(x) a |, \z(u) \, \y a \ < 1, Vfc, |a|, |/3| < 2r. 

From this, and with l(r) := 2Z(r), we immediately deduce that |z 7 | < 1 whenever 
|-y ] < 2r. In particular, L y (H) + L z (h) > — J^p \Hp\ ~ S 7 I^Ij which proves that 
inf Qj( r ) > — oo, and so inf Q r > — oo for r sufficiently large. 
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Let p := inf P = minP (by Theorem l2.3p . let r > l(ro), and let (z r ,y r ) be a nearly 
optimal solution of Q r with value 



(5.8) inf Q r < L z r(h) + L yr (H) < inf Q r + - l< p + 

Complete the finite vectors y r and z T with zeros to make them infinite sequences. 
Since for arbitrary s € N one has |y£|, \zZ\ < 1 whenever \a\, \j\ < 2s, provided r 
is sufficiently large, by a standard diagonal argument, there exists a subsequence 
{rfc} and two infinite sequences y and z, with |y Q | < 1 and |z 7 | < 1, for all a, 7, 
and such that 

(5.9) lim y r a » = y a Va£ N n ; lim 2^ = z 7 V7 € N X N™ X N m . 

— >oo k — >oc 

Next, let r be fixed arbitrarily. Observe that M rk (y Tk ) y implies M r (y rk ) y 
whenever rj, > r, and similarly, M r [z Tk ) y 0. Therefore, from (|5.9p and M r (y Tk ) > 
0, we deduce that M r (y) >z 0, and similarly M r (z) >: 0. Since this holds for 
arbitrary r, and \y a \, \&y\ < 1 for all a, 7, one infers from Proposition 13.21 that y 
and z are moment sequences of two measures v and p with support contained in 
[—1, 1]™ and [0, 1] x [—1, 1]" x [—1, l] m respectively. In addition, from the equalities 
yQ k = 1 and z^ k = 1 for every k, it follows that v and p are probability measures 
on [-1, 1]", and [0, 1] x [-1, 1]™ x [-1, l] m . 

Next, let (t, a) 6 N x N™ be fixed, arbitrarily. From 

V* - 5(0, x ) - L z r k {dg/dt + (V x g, /)) = 0, with g = (^x Q ), 

and the convergence (|5.9|) , we obtain 

L y ( ffl ) - 5 (0, i ) - ^(%/at + (V x g, /)) = 0, with g = {t*>x a ), 

that is, (Cg, {p,v)) = (g, <5(o,x ))- Since (t, a) G N x N™ is arbitrary, we have 

(g,£*(p,is)) = (Cg,(p,v)) = (g,5( 0}Xo )) Vg6l[f,4 

which implies that C*{p, v) = Sm xo y 

Let z(x), z(u) and z(t) denote the moment vectors of the marginals of p with 
respect to the variables x £ K", u S R m and t £ R, respectively, i.e., 

z{x) a = [ x a p(d{t,x,u)) Vaef, z{uf = [ u p p(d(t,x,u)) V/3 e N m , 



and z(<)' c = / i fc jit(cZ(t, sc, u)) for every k £ N. 

With r fixed arbitrarily, and using again (|5.9[) , we also have M r (9jy) y for 
every j 6 Jt , and 

M r (t!j z(x)) b VjeJ, M r (to t z(!i))H Vfc € W, M r (t(l - 1) z(t)) h 0. 



Since X, K and U satisfy Putinar's condition (see Definition [33]), from Theorem 
(Putinar's Positivstellensatz), y is the moment sequence of a probability measure 
v supported on K C [—1, 1]™. Similarly, z(x) is the moment sequence of a measure 
p x supported on X C [—1, 1]™, z(u) is the moment sequence of a measure p u 
supported on U C [— l,l] m , and z(t) is the moment sequence of a measure p 1 
supported on [0,1]. Since measures on compact sets are moment determinate, 
it follows that p x , p u , and p l are the marginals of p with respect to x, u and t 
respectively. Therefore, p has its support contained in S, and from C*(p, v) = 5(o iXo ) 
it follows that (p, v) satisfies all constraints of the problem P. 



24EAN B. LASSERRE, DIDIER HENRION, CHRISTOPHE PRIEUR, AND EMMANUEL TRELAT 



Moreover, one has 

lim infQ rfe = lim L z r k (h)+L y r k (H) (by Q) 

&j— >-oo k— too 

= L z (h) + L y (H) (byiD) 

= Jhdfi + J H dv </5 = minP. 

Hence, {fj,,v) is an optimal solution of P, and minQ, r j minP (the sequence is 
monotone nondecreasing) . Item (i) is proved. 
Item (ii) follows from Theorem 12.31 (hi) . □ 

5.6. Proof of Theorem 14.11 It suffices to prove that v p — > v as p — > +oo. For 
every integer p, v p = minP p is attained for a couple of measures (ji p , u p ). As in the 
proof of Theorem 12. 31 the sequence {(fi p , f p )} P 6N is bounded in A^(S) + x M(K) + , 
and hence, up to a subsequence, it converges to an element (/i, v) of this space for 
the weak * topology. 

On the one hand, by definition, C*(fX p , v p ) = S(p,x ) f° r every p. On the other, 
£* tends strongly to £*, and so C*(fi,v) — 5(o, Xo )- Moreover, since (h p ,H p ) tends 
strongly to (h,H) in Ci(S) x Ci(K), one has 

V P = ((Ppr v v)> ( h P> H p)) y V )' ( h > H ))> 

and so v < ((/i, v), (h, H)). We next prove that v = {(p, is), (h, H)). 
Since v p ) is an optimal solution of P p , 

{{v-p,v p ),(h p ,H p )) < ((p,,D),(h p ,H p ), V(/2,P) | C*(p,,u) = S^ Xo y 

Hence, passing to the limit, 

((/i, v), (h, H)) < ((fi,i?),(h,H), V(fi,v) \ C*(fl,v) = 8(o, Xo) , 

and so, (//, v) is an optimal solution of P, i.e., v = ((/i, v), (h, H)). □ 

Acknowledgments 

The research of Didier Henrion was partly supported by Project No. 102/06/0652 
of the Grant Agency of the Czech Republic and Research Program No. MSM6840770038 
of the Ministry of Education of the Czech Republic. 

References 

[1] E.J. Anderson, P. Nash. Linear Programming in Infinite-Dimensional Spaces. John Wiley, 
Chichester, UK, 1987. 

[2] R. Beals, B. Gaveau, P.C. Greiner. Hamilton-Jacobi theory and the heat kernel on Heisenberg 

groups. J. Math. Pures Appl. 79 (2000), 633-689. 
[3] C. Berg. The multidimensional moment problem and semigroups. Proc. Symp. Appl. Math. 

37 (1987), 110-124. 

[4] V. Borkar. Convex analytic methods in Markov decision processes. In: Hanbook of Markov 
Decision Processes. E.A. Feinberg and A. Shwartz, Eds. Kluwer Academic Publishers, 2002, 
377-408. 

[5] O. Bokanowski, S. Martin, R. Munos, H. Zidani. An anti-diffusive scheme for viability prob- 
lems. Submitted. 

[6] R.W. Brockett. Asymptotic stability and feedback stabilization. In: Differential geometric 
control theory. R.W. Brockett, R.S. Millman and H.J. Sussmann, Eds. Birkhauser, Boston, 
1983, pp. 181-191. 

[7] I. Capuzzo-Dolcetta, P.L. Lions. Hamilton-Jacobi equations with state constraints. Trans. 
Amer. Math. Soc. 318 (1990), 643-683. 



OPTIMAL CONTROL AND LMI-RELAXATIONS 



25 



[8] Cho, Moon Jung, R.H. Stockbridge. Linear programming formulation for optimal stopping 
problems. SIAM J. Control Optim. 40 (2002), 1965-1982. 

[9] C. Coatmelec. Approximation et interpolation dcs fonctions diffcrentiablcs de plusieurs vari- 
ables. Annates Scientifiques de I'Ecole Normale Superieure, 3eme serie, 83, 4 (1966), 271—341. 
[10] D.A. Dawson. Qualitative behavior of gcostochastic systems. Stock. Proc. Appl. 10 (1980), 
1-31. 

[11] W.H. Fleming, D. Vermes. Convex duality approach to the optimal control of diffusions. 
SI AM J. Control Optim. 27 (1989), 1136-1155. 

[12] R. Fletcher. Practical methods of optimization. Vol. 1. Unconstrained optimization. John 
Wiley & Sons, Ltd., Chichester, 1980. 

[13] P.E. Gill, W. Murray, M.H. Wright. Practical optimization. Academic Press, Inc., London- 
New York, 1981. 

[14] R.F. Hartl, S.P. Sethi, R.G. Vickson. A survey of the maximum principles for optimal control 
problems with state constraints. SIAM Rev. 37 (1995), 181-218. 

[15] K. Helmes, R.H. Stockbridge. Numerical comparison of controls and verification of optimality 
for stochastic control problems. J. Optim. Theory Appl. 106 (2000), 107-127 

[16] K. Helmes, S. Rohl, R.H. Stockbridge. Computing moments of the exit time distribution for 
Markov processes by linear programming. Oper. Res. 49 (2001), 516-530. 

[17] D. Hcnrion, J.B. Lasserre. Solving nonconvex optimization problems. IEEE Control Systems 
Magazine 24, 2004, pp. 72-83. 

[18] O. Hernandez- Lerma, J.B. Lasserre. Discrete- Time Markov Control Processes: Basic Opti- 
mality Criteria. Springer Verlag, New York, 1996. 

[19] O. Hernandez-Lerma, J.B. Lasserre. Further Topics in Discrete- Time Markov Control Pro- 
cesses. Springer Verlag, New York, 1999. 

[20] O. Hernandez-Lerma, J.B. Lasserre. Markov Chains and Invariant Probabilities. Birkhauscr 
Verlag, Basel, 2003. 

[21] O. Hernandez-Lerma, J.B. Lasserre. Approximation schemes for infinite linear programs, 
SIAM J. Optim. 8, 1998, pp. 973-988. 

[22] O. Hcrnandcz-Lcrma, J.B. Lasserre. The linear programming approach, In: Hanbook of 
Markov Decision Processes. E.A. Feinbcrg and A. Shwartz, Eds. Kluwer Academic Pub- 
lishers, 2002, pp. 377-408. 

[23] D. Hernandez-Hernandez, O. Hernandcz-Lcrma, M. Taksar. The linear programming ap- 
proach to deterministic optimal control problems. Appl. Math. 24, 1996, pp. 17-33. 

[24] M.W. Hirsch. Differential topology. Graduate Texts in Mathematics 33, Springer- Verlag, 
1976. 

[25] D. Jacobson, M. Lele, J. L. Speyer. New necessary conditions of optimality for control prob- 
lems with state-variable inequality constraints. J. Math. Analysis Appl. 35 (1971) 255-284. 
[26] J.L. Krivine. Anneaux preordonnes. J. Anal. Math. 12 (1964), 307-326. 

[27] T.G. Kurtz, R.H. Stockbridge. Existence of Markov controls and characterization of optimal 

Markov controls. SIAM J. Control Optim. 36 (1998), 609-653. 
[28] J.B. Lasserre. Global optimization with polynomials and the problem of moments. SIAM J. 

Optim. 11, 2001, pp. 796-817. 
[29] J.B. Lasserre, T. Prieto-Rumeau. SDP vs. LP relaxations for the moment approach in some 

performance evaluation problems. Stoch. Models 20 (2004), 439-456. 
[30] J.B. Lasserre, T. Prieto-Rumeau, M. Zervos. Pricing a class of exotic options via moments 

and SDP relaxations. Math. Finance 16 (2006), 469-494. 
[31] J.B. Lasserre, C. Prieur, D. Henrion. Nonlinear optimal control: Numerical approximations 

via moments and LMI relaxations. Proc. 44th IEEE Conference on Decision and Control, 

Scvilla, Spain, December 2005, pp. 1648-1653. 
[32] H. Maurer. On optimal control problems with bounded state variables and control appearing 

linearly. SIAM J. Cont. Optim. 15 3 (1977) 345-362. 
[33] J. Nash. The imbedding problem for Ricmannian manifolds. Ann. Math. 63 (1956), 20-63. 
[34] H.J. Pesch. A practical guide to the solution of real-life optimal control problems. Control 

Cybernet. 23, 1-2 (1994), 7-60. 
[35] L.S. Pontryagin, V.G. Boltyanskij, R.V. Gamkrelidze, E.F. Mishchcnko. The Mathematical 

Theory of Optimal Processes John Wiley & Sons, New York, 1962.. 
[36] C. Prieur, E. Trelat. Robust optimal stabilization of the Brockett integrator via a hybrid 

feedback. Math. Control Signals Systems 17, 3 (2005), 201-216. 



28EAN B. LASSERRE, DIDIER HENRION, CHRISTOPHE PRIEUR, AND EMMANUEL TRELAT 



[37] M. Putinar. Positive polynomials on compact semi-algebraic sets. Ind. Univ. Math. J. 42, 
1993, pp. 969-984. 

[38] K. Schmiidgen. The X-moment problem for compact semi-algebraic sets. Math. Ann. 289 
(1991), 203-206. 

[39] M.H. Soner. Optimal control with state-space constraints. SIAM J. Control Optim. 24 (1986), 
552-561. 

[40] J. Stocr, R. Bulirsch. Introduction to Numerical Analysis. Third Edition, Springer- Vcrlag, 
New York, 2002. 

[41] O. von Stryk, R. Bulirsch. Direct and indirect methods for trajectory optimization. Ann. 
Oper. Res. 37, 1992, pp. 357-373. 

[42] E. Trelat. Controle optimal: theorie et applications. Vuibert, Paris, 2005. (in french) 

[43] F.-H. Vasilescu. Spectral measures and moment problems. Spectral Theory and Its Applica- 
tions, Theta 2003, 173-215. 

[44] R. Vinter. Convex duality and nonlinear optimal control. SIAM J. Control Optim. 31, 2 
(1993), 518-538. 

LAAS-CNRS and Institute of Mathematics, University of Toulouse, France. 
E-mail address: lasserreSlaas.fr 

LAAS-CNRS, University of Toulouse, France and Czech Technical University in 
Prague, Czech Republic. 

E-mail address: henrlonaiaas.fr 

LAAS-CNRS, University of Toulouse, France. 
E-mail address: cprieuraiaas.fr 

University of Orleans, France 

E-mail address: emmanuel.trelatauniv-orleans.fr 



